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ABSTRACT 

This paper investigates the spectrum of the iteration operator of some finite element pre- 
conditioned Fourier collocation schemes. The first part of the paper analyses one-dimensional 
elliptic and hyperbolic model problems and the advection-diffusion equation. Analytical ex- 
pressions of the eigenvalues are obtained with use of symbolic computation. The second 
part of the paper considers the set of one-dimensional differential equations resulting from 
Fourier analysis (in the tranverse direction) of the 2-D Stokes problem. All results agree with 
previous conclusions on the numerical efficiency of finite element preconditioning schemes. 
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Science and Engineering (ICASE), NASA Langley Research Center, Hampton, VA 23665. 
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1. INTRODUCTION 


In the recent past, Chebyshev collocation schemes have been applied extensively to the 
numerical integration of the Navier-Stokes equations [1, 3, 4]. For scalar elliptic problems, it 
is well known that the condition number of the matrix system of discrete algebraic equations 
increases rapidly with N, the number of degrees of freedom of the problem at hand. Therefore, 
the preconditioning technique seems to be the only adequate tool in order to overcome 
this numerical burden. The present authors [5, 6] demonstrated that finite elements (FE) 
constitute powerful preconditioners for general second-order elliptic equations. In [3], several 
fluid flow elements in velocity-pressure formulation were investigated. From the analysis of 
the eigenspectrum of the iteration operator, it was shown that the Q2-Q1 element is the 
best choice for the steady Stokes problem. As all the previous analyses on finite element 
preconditioning were carried out numerically, the present note aims at analytical results 
through use of symbolic manipulation languages (cfr. [10]). 

For the sake of simplicity, we will restrict ourselves to the study of a finite element 
preconditioned Fourier collocation scheme. In this case, the mesh size is uniform and the 
algebra is considerably reduced. In Section 2, a one-dimensional model is considered. The 
collocation process is preconditioned by Lagrangian linear, quadratic, cubic and Hermite 
cubic elements, respectively. The Richardson iteration method is set up with these FE 
preconditioners as approximate operators and algebraic solvers. Using the spatial structure 
of the eigenvectors of the Fourier solutions, one may perform a full analysis of the eigenvalues 
of the iteration operator. This theoretical investigation corroborates the previous numerical 
analyses [6]. In Section 3, a one-dimensional hyperbolic model is investigated using linear and 
quadratic FE preconditioning. The upwinding technique is also examined. A further model 
consists in an advection-diffusion equation. In Section 4, the Stokes problem is reduced to 
a 1-D incompressible flow model amenable to Fourier discretization. The Q2-Q1 and Q1-P0 
elements are candidates as preconditioners. A similar Fourier analysis is done. The results 
corroborate numerical experiments carried out in the framework of Chebyshev collocation 


[3]. 


2. ELLIPTIC MODEL 

Let us first consider the simple elliptic problem: 

-u xx = /, 0 < x < 2ir, (2.1) 

with periodic boundary conditions. The subscript indicates partial derivative. The Fourier 
approximation of the dependent variable u is 

N/ 2-1 

U N = ± V ipx, > 0<j<N, (2.2) 

P=- N/2 

where Up are the discrete Fourier coefficients and Xj the collocation points defined by 

X, = 22 , j € [O.JV[. (2.3) 
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The linear system corresponding to (2.1) may be found in [1] and will be denoted by L c . 
The eigenfunctions of (2.1) are 

tj(p) = 0 <j<N, (2.4) 

with the corresponding eigenvalues 

A(p)=p 3 , p€[-y,y-l]. (2.5) 

The collocation problem will be preconditioned by finite elements. Introducing the approxi- 
mate FE operator L, the preconditioned Richardson iteration is written as: 

u k+1 = u k - ot k L~ l (L c u k - /), (2.6) 

where k is an iteration index, a k a relaxation factor and u, f the vectors corresponding to the 
unknowns and source terms at the collocation points. The convergence of (2.6) is governed 
by the spectral radius p(A) of the iteration operator defined by A = I — a k L~ 1 L c . The 
optimal value of the relaxation factor is: 


&opt 


+ A„ 


(2.7) 


where A mtn and A moi are the minimum and maximum eigenvalues of L~ l L c . An approximate 
estimate of the number of iterations n needed to reduce the error norm by a factor £ is given 
by 

n = - log £/£«,( A), (2.8) 

where Roo(A) = — log p{A) is the asymptotic rate of convergence of the iteration matrix. 
The spectral radius p{A) which is involved in the error reduction process with the use of 
(Eq. (2.7)) is given by 

p(A) = (2.9) 

^ max i 

In order to investigate this quantity for various preconditioners, we have to define the 
finite element problem more precisely. 

Lagrangian linear elements, Hermite cubic elements (i.e., Ql } P 3 in Ciarlet’s notations [2]) 
as well as higher-order Lagrangian interpolants have their vertices at the Fourier collocation 
grid (2.3). However, for Lagrangian quadratics ( Q2 ), mid-points are added at 

2tt 1 

= Jfti + 2 } ’ j e [0 ' N[ ' (2 - 10) 

while for Lagrangian cubics (Q3), we have additional grid nodes located at 


x j+ ^,x j+ i, i€[0,JV[, (2.11) 

with obvious definitions. For the Lagrangian case, the FE unknowns are the nodal values, 
while for the Hermite case, the unknowns are Uj and u'-, the prime denoting the first-order 
derivative of u. In [5], the iteration operator is written as: 

A = I- S^MnhLc, (2.12) 
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where S h is the stiffness matrix, M h the mass matrix, I h an interpolation matrix evaluating 
the Fourier interpolant of the collocation operator at the FE nodes. The mesh size of the 
FE grid is defined by 

h = 2 -k/N. 

For Q\ elements, h reduces to the unit matrix; for higher-order interpolants however, the 
structure of this matrix is more complicated. In order to avoid writing the details of Ih, 
we will systematically assume the use of static condensation. Consequently, the iteration 
operator may be written: 

A — I — §h'M h L c , (2.13) 

where §h and Mh refer to stiffness and mass matrices after static condensation. 


2.1. Linear Lagrangian Elements 

For an interior node, the expressions for the stiffness and mass matrices are well-known: 


ShUj = \ ( ' ~ u i~ i + 2u,- - u i+1 ) , (2.14) 

M h fj = ^(/;-i + 4 fj + (2.15) 

Fourier analysis of (2.14), (2.15) with the eigenfunction (2.4) leads to the expression of the 
eigenspectrum of S^MhLc, denoted by a (p), 


tf(p) = 


(ph/2) 3 (2 + cos ph) 
sin 3 ph/2 3 


-£<,<£-i. 
2 ~ 2 


(2.16) 


Typically, the second factor in the right-hand side of (2.16) comes from the contribution 
of the mass matrix. In the case of finite difference (FD) preconditioning, this factor is 
one. For p = 0, a(p) = 1, while for p = -N/2 ) a{p) = tt 2 /12. This last value should be 
compared to the FD equivalent which is c{p) = 7t 3 /4 [7]. The eigenvalue spectrum of the 
FE preconditioning is reduced because of the beneficial presence of the mass matrix. Fig. 1 
shows the behavior of cr(p) with respect to p for h = 27r/100. The function has a minimum 
value equal to 0.693. Therefore, the optimum value for os is 


&opt — 1*18, 


(2.17) 


and over-relaxation is possible for FE preconditioning unlike the FD preconditioning where 
under-relaxation is required to converge. In practice, the Q1 preconditioning with a spectral 
radius of 0.18 converges twice as fast as the FD preconditioner whose spectral radius is of 
the order of 0.42. 
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2.2. Quadratic Lagrangian Elements 

The equations related to nodes j and j i | may be cast in the following matrix form: 


0 


_8 16 -8 

3/ * * 0 0 -8 16 -8 


-8 14 -8 


0 

1 


Uj - 1 N 


, j 

(fi- 1\ 

u i~h 

h I 

\ 8 1 0 M 

fj~\ 

Uj 

~ 15 ( 


fi 

Uj+l J 

S' 

O 

O 

1—4 

OO 
1— * 

fi + i J 



The use of static condensation eliminates the contribution of u 4 
reduces to only one relationship for node j on the collocation gri 


(2.18) 

_i and u j+ i and Eq. (2.18) 
ria 

(2.19) 


1 + 2u i ~ u i+i) = + fi + fj+ i)- 

Let us notice that for <3/i on the left hand-side of (2.19), we recover the stiffness matrix Sh, 
associated to Q1 elements whereas in the right-hand side, Mh corresponds to a different 
quadrature rule. Carrying out the Fourier analysis of (2.19), one obtains 


“ sin 2 ph/2 3 ^ + 2 cos W 2 ))’ "J - p - ~2 ~ L 


( 2 . 20 ) 


For the particular values p = 0 and p = -N/2, a(p) is equal to 1 and 7r 2 /12,respectively. As 
a{p) is a monotically decreasing function with respect to p (Fig. 1), the optimum value of a 
is 

«opt = 2/(1 + 7T 3 /12) = 1.0974, (2.21) 

and the corresponding spectral radius p(A) is equal to 0.0974. 

2.3. Cubic Lagrangian Elements 

For the sake of compactness, we give the local stiffness and mass matrices over the uniform 
mesh: 


( 2 . 22 ) 


/ 

r 37/20 

-189/80 

27/40 

-13/80 \ 

2 

-189/80 

27/5 

-297/80 

27/40 


h 

27/J40 

-297/80 

27/5 

-189/80 

\ 

k -13/80 

27/40 

-189/80 

37/20 

J 


16/105 

33/260 

-3/70 

19/840 \ 


h 

33/280 

27/35 

-27/280 

-3/70 


: 2 

-3/J70 

-27/280 

27/35 

33/280 



19/840 

-3/70 

33/280 

16/105 j 



(2.23) 


Assembling the matrices of (2.22), (2.23) over two adjacent elements and eliminating the 
four unknowns attached to the interior nodes Uy±i/3,u,±2/ 3 , we are left with the relation: 

+ 2 Uj - u i+a ) = ^(^ + 4/,-f + + 3 4+i + \fj + f + ^r 1 )- (2-24) 
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Here again, as in the previous case, we recover in the left-hand side of (2.24) the stiffness 
matrix of Ql elements. Fourier analysis of (2.24) leads to the eigenvalue spectrum cr(p) : 


(ph/2) 2 1,4 jp/i ph 2 N N , . 

"W = Bin* (pft/2 ) 10 J3 C ° S T + 3c ° S (2 ’ 25) 


The particular values of a(p ) corresponding to p = 0 and p = —NJ 2 are cr(p) = 1 and 
497 t 2 /480 ~ 1.007522. The optimal value, of the relaxation factor is given by 2/(1 + 
497r a /480) ~ 0.9963 and the corresponding spectral radius p(A) is equal to 0.00375. 

Looking back at the results in the previous subsections, one observes that the spectral 
radius p(A) diminishes with increasing polynomial degrees. This does not mean however 
that one should use higher order elements in the preconditioning because they involve more 
computational work as the bandwidth of the algebraic system increases. 


2.4. Cubic Hermite Elements 

At node j, the discrete equations are 

-^K- 1 - 2u > + u ;+0 - t^K-i - u i+i) = - /?+i) 
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420 

h 9 26 9 

+ 7 { lQ fj - 1 + ~5 fi+ lQ fj+l) ' 


(2.26) 


- Wj+i) - 1 - 8 u 'i + u j+i) 


Fourier analyzing (2.26), one gets the spectrum: 


13 


420 

h 3 

140 


A’C/i-x - /,«) 

- \n + /;«)■ 


(2.27) 


(ph/2) 


13 


= 6 .i n\ P h/2) - ( P A/2) sin(pft/a) c«&W 7 (26 + ,C “* + 

|p/i|e[0,4 (2.28) 

For p = 0 and -N/ 2, <r(p) = 1 an 177 t 2 /168 ~ 0.9987, respectively. Fig. 1 displays the 
behavior of a(p), which is first slightly decreasing with respect to p achieving its minimum 
value 0.97722 and then increasing to 1. The optimal value of a is 1.01152 and the cor- 
responding spectral radius 0.0115, an intermediate value in between those of Q2 and Q3 
preconditioning. 

3. HYPERBOLIC PROBLEMS 

We now turn our attention to the first-order differential equation 

u x = f, (3.1) 
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in the periodic case. The eigenfunctions of (3.1) are again (2.4) with the eigenvalues 

, . , N N ,, 

A(p) = tp, P€[— .y-lj. 


3.1. Linear Lagrangian Elements 

Using the standard Galerkin approach, a centered scheme is produced and yields the 
discrete equation 

— ~ U, '~* = £(/i-. + ifi + /;+.)• (3.2) 

By Fourier analysis, we obtain: 

4)4^, Me [ 0,4 (3.3) 

sinphf o 

For p = 0, a(p) = 1, while cr(—N/ 2) is obviously unbounded as in the FD case. The presence 
of the mass matrix does not help to circumvent the difficulty. 

One may proceed, however, using an upwinding technique. This has been a key step 
to treat hyperbolic problems. In finite elements, the method useB separate test and trial 
function spaces, i.e., the Petrov-Galerkin method. There are several ways to implement 
upwinding. Let us introduce the weight functions iy’(r),(i = 1,2) defined on the reference 
interval [— 1 , 4*1] by Heinrich and Zienkiewicz [8]: 

w\r) = <p*'(r) + (-l)’eF(r), -1 < r < 1, (3.4) 


where <p*(r) are the linear Lagrangian trial functions, F(r) an auxiliary quadratic element 
vanishing at both end points 

F(r) = j(l - r’), (3.5) 

and e the upwinding parameter to be given independently. Using <p l and w' as trial and test 
functions respectively, one gets — 


1 + e e 1 — e 

~2 + i Ui + - w ** 1 


2 -f- 3e 
12 


fj-i + g/i + 


2 — 3e 
12 ^ +1 ‘ 


(3.6) 


This equation reduces to (3.2) when e = 0 (no upwinding). With a value of e as yet undefined, 
the Fourier analysis of (3.6) leads to complex eigenvalues given by: 

o-(p) = r-: -(3ephsmph + 2iph(2 + cos ph)), \ph\ 6 [0,7r]. (3.7) 

6 e(l — cosp/i) + isinp/i 


For p = 0, cr(p) = 1 while for p = — N/2,a(p ) = i/2 independently of the value of e. 
Upwinding forces the spectrum of S^MhLc almost entirely inside the unit circle, as shown 
on Fig. 2 where the eigenvalues (3.7) have b een plo tted for e = 1. The s pectral radius 
of the matrix A in this case is equal to \/5/2 and under-relaxation is required in order to 
ensure convergence of the preconditioning iterations. The eigenvalues (3.7) being complex, 
the evaluations of a ^ and the spectral radius p(A) are no longer given by (2.7) and (2.9). 
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3.2. Quadratic Lagrangian Elements 


Another way to solve (3.1) consists in using a FE preconditioner based on quadratic 
Lagrangian elements. Applying the Galerkin approach and assembling the contributions of 
two adjacent elements at node j, one obtains a set of three equations related to nodes j and 
j ± 5 similar to (2.18), which are cast in matrix form: 


-I 0 


2 

3 

0 


0 


0 

2 

3 

0 


vj- 1 \ 






( fj-l ^ 


h 

2 

16 

2 

0 0 \ 


Uj 

= 30 

_1 

2 

8 

2 _1 

fi 

Uj + 1 / 

l o 

0 

2 

16 2 J 

/; + i > 


(3.8) 


If static condensation is carried out through Gaussian elimination, the contribution of the 
exterior nodes Vj-t,Vj+t disappears and one is left with a staggered scheme: 


|(u J -_i - u j+ 1) = ^(-/j-i + + 18 /i + 12 /j+$ “ fi+ 0> 


S' 60 

for Fourier analysis. Its spectrum is given by: 

9 + 12 cos ^ — cos ph 


a(p) = 


Eh 

2 


sin y 


20 


\ph\ G [0 ,tt]. 


(3.9) 


(3.10) 


It is monotonically decreasing and bounded by cr(0) = 1 and a(—N/ 2) = 7r/4 as shown on 
Fig. 3. One should notice that the first term in the right-hand side of (3.10) is identical 
to the spectrum obtained in the FD case where the function is computed on the main grid 
while the derivative is evaluated by first-order differences on a staggered grid. The second 
term whose maximum value is equal to 1 is induced by the presence of the mass matrix and 
reduces to unity in the FD case. The optimal value of a is equal to 


a opt = 2/(1 + 7r/4) = 1.1202, 

and the corresponding spectral radius p(A) is equal to 0.1202. This staggered scheme gener- 
ated by Q2 elements is the key of success for FE preconditioning of Navier-Stokes problems. 
This excellent behavior explains the reason why in Demaret-Deville [5], the relaxation pa- 
rameter was almost independent of the Reynolds number. 


3.3. Advection-Diffusion Model 

The last scalar model analyzed in this paper is the one-dimensional advection-diffusion 
problem. The differential equation writes 

— ku xx + cu x = f(x), (3-11) 

where k is the diffusion coefficient and c the constant advection velocity. Particular interest 
bears on advection dominated problems which impose severe conditions on the element mesh 
size (cfr. Thomasset [9]). With the eigenvectors (2.4), the eigenvalues of (3.11) are 

N N 


where 7 is the cell Reynolds number defined by 7 = chj k. 

Using the linear FE basis and upwinding introduced in the previous section, one gets the 
discrete equations 


-(1 + + (2 + - (1 - 7—5— ) u i+i = Yjl( 2 + + 8/, + (2 - 3e)/,+,|. 

( 3 - 12 ) 

where e is the upwinding parameter of Eq. (3.4). With no upwinding (i.e., e = 0), stability 
requirements restrict 7 to values < 2. The Fourier analysis of (3.12) is straightforward. The 
eigenvalues of (3.12) are complex and given by: 


a (p) = 


^(2 + cos ph) + sin ph + i[7^£(2 + cos ph) — e^^-sinp/i] 


2(2 + ey) sin 2 ^ + iy sin ph 


, \ph\ e [o,tt]. 

(3.13) 


In absence of upwinding, Eq. (3.13) reduces to an analytical expression whose real and 
imaginary parts may be written in compact form: 

n . . .. , 4p/i sin 2 ^ + 7 2 sinph 2 + cos ph 

Re(a{p)) = ph — . f — — 7 5 , 

16sm +7 a sinpn 0 


Im(a(p)) = 7 ph 


4 sin 2 *y — ph sin ph 2 + cos ph 
16 sin 4 y + 7 2 sinph 3 


(3.14) 


The factor (2 + cospA)/3 in the right-hand side of (3.14) is another example of the contri- 
bution of the mass matrix in FE preconditioning. Like in Eq. (2.16), this factor reduces 
to unity in the expression of the eigenvalues corresponding to FD preconditioning. One 
can draw similar conclusions to the diffusion problem, except for the complex nature of the 
eigenvalues. Introducing ph = 0 and ph = it into the eigenvalues of the FD case gives the 
bounds of the spectrum: 


1 < Re(a FD ) < 0 < Im{a FD ) < ^-. (3.15) 

In the FE case, the upper bounds are reduced by a factor 3 because of the presence of the 
mass matrix. Figure 4 displays the result (3.14) for both FD and FE preconditionings and 
for two different values of 7 (i.e., 0.2 and 2). Even in the limit case 7 = 2, the spectrum of 
A for FE preconditioning lies inside the unit circle. Reducing the value of the cell Reynolds 

number brings the eigenvalues closer to the real axis. 

Figure 5 exhibits (3.13) with 7 = 2, with and without upwinding. In the upwinding 
case, e was chosen equal to 1. The eigenspectrum is rotated counterclockwise and slightly 
stretched inducing a somewhat larger spectral radius. 


4. STOKES EQUATIONS 

Let us write the Stokes equations in stress formulation: 

diva. + p£_ = 0, 


(4.1) 
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divv = 0. (4-2) 

The symbol g. denotes the stress tensor/ p is the density, / the body force term and v is 
the velocity field. Eq. (4.1) is the momentum equation and Eq. (4.2) enforces the continuity 
constraint. The 2-D Stokes problem may be reduced to a 1-D problem if the solution of 
(4.1), (4.2) is sought as a Fourier mode: 

■H = u{x)e ikv , p = p(x)e ikv . (4.3) 


Introducing (4.3) in (4.1), (4.2), we get: 


dv du 2 ... dv . . . 

+ 2 "a? + ,kK ' ku + ai ) + pU = 0 ’ 

d dv 

p-^{ik u + -^) - ikp - 2pk 2 v + pf y = 0, 


du .. . 

— +ikv = 0, 
Ox 


0 < x < 2r. 


(4.4) 

(4.5) 

(4.6) 


The velocity and pressure fields are assumed to be 27r-periodic. This 1-D problem is dis- 
cretized in the x direction using Fourier series of type (2.2) for each variable. The discrete 
collocation equations are preconditioned by finite elements such as the Q2 — Q 1 and Ql — P 0 
elements. The FE equations come from Galerkin projection. Introducing ipi and <pt the trial 
functions for the FE approximations of the velocity and pressure fields, respectively, such 
that 

M N 

v(x) = Y^*> p( x ) = ( 4 - 7 ) 

i=i i=i 

the discrete FE equations are obtained by use of the divergence theorem as tool for the 
integration by parts with the notation /* = /, f y = g- 


X)[2 pAji + pk 2 Bj{\ui - ikp Y C H v i ~ 12 D i& 1 = 12 Q <j < M, 

i ill 

ikp Y Cji u l + 2pk 2 Bp + pAji]v t + ikY E flPi = 12 B ji9h 0 < j < M, 

l l il 

-Y D l u < ~ ^Y E h = °» 0<j<N. 

l l 

In (4.8-4.10), the various matrices are defined by the relationships: 


A * = JltT& dx ’ = c >‘ = 


din 


Dp = J Eji = J (pjipidx. 


(4.8) 

(4.9) 
(4.10) 


(4.11) 
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4.1. Q2-Q1 Elements 


For this element, M = 2 N. Carrying through the algebra involved by the quadratures 
(4.11) and assembling by direct stiffness the contributions of the two elements connected to 
node j, we obtain: 


2a f 

j-i _ Suy-i/a + 14 Uj — 8uj + i/ 3 + uy+i) 

+ + 8u > + 2u J+i/a ~ u i+i) “ 1 “ 4vy-i/a + 4u J+1/a - v i+1 ) 


-g(Pi-i - Pj+ i) = 2Q (-/j-i + 2/j-i/a + 8 fj + 2fj + i/ 2 - fj+i), 

iu iik 2 h._ 2 ikii , . 

— (-4uj + 8uj + i/a - 4u J+ i) + -jjQ-(2uy + 16ii,+i/3 + 2u J+1 ) — (-Vj + uy+i) 

“g(Pi - Pj+ 1) = gQ( 2/i + 16/,+i/a + 2 /y+i), 


(4.12) 

(4.13) 


ikfi. . ak 2 h, „ _ . 

— C — "“i— i + 4u J -_i/a — 4 u j+ i/ 2 + “j+i) jg”(— ty-i + 2uj_i/a + 8 Vj + 2v J+1 / 2 — uy+i) 


2a, ikh 

+y( v i-i - 8uj_i/a + 14v y - 8v,- +1/a + vy+i) - —p, 

h 

= 2q( — 5j-i + 2^,-i/a + 8^j + 2g J+1 / a - <7y+i), 

(4.14) 

2ika . . 8u. 2uk 2 h. „ . 

-jj-W - u i+i) + 3^(- u i + 2u i+i/3 - u i+i) 15 - + 8 v 3'+i/ 2 + v i+i) 

(4.15) 

~~^~(Pj + Pi+i) = J^(f) + ^fj+}/2 + fj+i), 

1 2ik 

-(— — 4u J _i/a + 4u, + i/a + u J+ i) + -^-(uj-i/a + vy + Vj+ l/a) = 0. (4-16) 

Eqs. (4.12), (4.13) and (4.16) correspond to momentum and incompressibility relations, 
while (4.13) and (4.15) are the momentum equations associated to mid-node ccj+j/a. Similar 

expressions hold for mid-node sy_i/a with appropriate shifts for the indices. 

Now, static condensation represents a formidable task and is greatly helped by the sym- 
bolic manipulation program. Elimination of «y±i/a, uy±i/a leads to a matrix system of order 
3. 

The full Fourier solution gives the collocation matrix L c : 

( —2 fil 2 — k 2 fi —fikl —il \ 

L c = —fikl —fit 1 — 2fik 2 —ik I , (4-17) 

\ — il —ik 0 / 
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where l is the wavenumber in the x direction. 

The analytical computation of MhL c is performed as far as the symbolic program 
r an handle tractable expressions. Then numerical evaluation of the eigenspectrum is done. 
Because of the divergence-free constraint, a zero eigenvalue is systematically obtained. In 
Figures 6 and 7, the eigenspectrum of MhL c are plotted for two cases k = 1 and k = 
10, respectively. In these two figures, the lower curve shows the same behavior as the 
eigenspectrum of the elliptic problem solved by Q2 elements. For l = —Nf 2, cr(l) is equal to 
7r 2 /12. For the other curve, <j( 0) = 1 and a{—Nj 2) is close to 2.07. Therefore, the optimal 

a value is given by 

aopt = 2/(2. 07 + tt 2 /12) ~ 0.69, 


a value close to 2/3 obtained by Demaret-Deville [3] for a 2-D Chebyshev collocation dis- 
cretization of the Stokes problem preconditioned by Q2-Q1 elements. 


4.2. Q1-P0 Element 

The quadratures (4.11) provide less complicated discrete equations in this case: 

i -I- 2 Uj - Uj+i ) - + 4uj + u;+i) + 1 “ 

h 0 * ( 4 . 18 ) 

~Pj + Pj+i = g (/j-i + 4 /j + fj+ 1)> 

^p(u,_i - Uj + 1 ) - + 4 Vj + v,-+i) + t(~ v 3- 1 + 2v i - v i+0 

2 13 (4.19) 

+ Pi + 0 = + 49 J ' + 9j+1 ^ 

— + Uj + 1 + -^-(uj-i + 2vj + u J+ i) = 0. (4-20) 

Obviously, this element generates second-order differences for partial derivatives. When the 
mass matrix is involved, the standard weighted mean between three adjacent nodes appears 
in the expressions. No static condensation is needed in this case. Fourier analyzing Eqs. 
(4.18)-(4.20), the stiffness and mass matrices are now: 

/8gsin 3 (^) + ^x ii (2 + cos/i0 pk sin hi 2tsin(^) \ 

Sh = I pk am hi 4j| sin a (y) + ^-^(2 + cos hi) ikh cos(-j) I i 

\ 2% sin hi + cos W) 0 / 

M h = diag(^(2 + cos hi), ^(2 + cos hi), 0). (4.21) 

3 3 

In Figures 8 and 9, the eigenspectrum of S^ 1 MhL c are displayed for k = 1 and 10, respec- 
tively. In these two figures, the top curve is that of the elliptic model preconditioned by Q1 
element. The bottom curve starts from 1 for l = 0, decreases till a minimum value close to 
0.49 and then increases to reach a(—N/2) = 0.5. The optimum value is 

ocopt = 2/(1 + 0.5) = 4/3. 
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5. CONCLUSIONS 


In this paper, we have Fourier analyzed the eigenspectrum of the iteration operator for 
finite element preconditioning of Fourier collocation applied to one-dimensional problems. 
For elliptic models, this theoretical analysis confirms previous numerical findings, especially 
the beneficial presence of the mass matrix which reduces the bounds of the eigenspectrum. 
For first-order problems, linear elements without and with upwinding are considered. With 
quadratic elements, a staggered scheme is produced. Its eigenspectrum is bounded and ranges 
between 1 and 7r/4. Finally, a Stokes problem is reduced to a one-dimensional approach. 
Two types of elements are examined. The Q2-Q1 element leads to an optimum value of 
the relaxation parameter close to the value obtained by numerical analysis of preconditioned 
Chebyshev collocation. For Q1-P0 element, the method can be over-relaxed. 
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Figure 1. Eigenspectrum of S h x M h L c for an elliptic model. Preconditioners: Ql:*; Q2:*; 
Q3:A; Cubic Hermite: Q- 
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Figure 4. Eigenspectrum of an advection-diffusion problem. The curves outside the unit 
circle are concerned with FD preconditioning while the inside curves are related to FE 
preconditioning. Both top curves are obtained for the cell Reynolds number 7 = 2. The 
bottom curves are gotten for 7 = 0.2. 
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Figure 5. Eigenspectrum of the advection- diffusion problem preconditioned by Ql elements 
for 7 = 2, with (Q) and without (*) upwinding. The upwinding parameter is e = 1. 
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Figure 6. Eigenspectrum of the Stokes problem preconditioned by Q2-Q1 element. Here, 
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Figure 9. Same caption as Figure 8, but with k — 10. 
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